R Coding Workshop: 9th Meeting

GIS & Geospatial Data Analysis (Fall 2025)

GIS Coding Workshop: 9th Meeting

Outline

  • Data Import
  • ggplot

Data Import

download.file() download files from an url

grandison_url <- "https://www.gutenberg.org/cache/epub/13884/pg13884.txt"

download.file(url = grandison_url, destfile = 'grandison.txt', quiet = TRUE)
grandison_text <- readLines('grandison.txt')
grandison_text |> head()
[1] "The Project Gutenberg eBook of The History of Sir Charles Grandison, Volume 4 (of 7)"
[2] "    "                                                                                
[3] "This ebook is for the use of anyone anywhere in the United States and"               
[4] "most other parts of the world at no cost and with almost no restrictions"            
[5] "whatsoever. You may copy it, give it away or re-use it under the terms"              
[6] "of the Project Gutenberg License included with this ebook or online"

Components of ggplot

We can define an ggplot object by specifying:

  • data
  • aesthetic mapping
  • geometry layer

Patterns - use + to introduce more layers or components

bk_df <- read_csv('data/us-main-books.csv')
bk_df |> glimpse()
Rows: 48
Columns: 25
$ state                 <chr> "WY", "PA", "OH", "NM", "MD", "RI", "OR", "WI", …
$ division_code         <dbl> 8, 2, 3, 8, 5, 1, 9, 3, 4, 8, 5, 2, 7, 4, 4, 8, …
$ division_str          <chr> "Mountain", "Middle Atlantic", "East North Centr…
$ political_value_index <dbl> -19.7, 2.0, -0.7, 2.4, 8.5, 11.2, 4.0, 2.4, -10.…
$ median_income         <dbl> 4081.5, 4218.0, 2469.0, -2401.5, 19404.5, 8141.0…
$ hs_grad_rate          <dbl> 8.3380421, 2.3380421, 3.4380421, -0.6619579, 4.2…
$ college_grad_rate     <dbl> -2.1237299, -1.6237299, -2.9237299, -0.5237299, …
$ challenge_count       <dbl> 4, 148, 29, 3, 5, 3, 118, 10, 6, 0, 13, 25, 5, 1…
$ explicit_count        <dbl> 1, 15, 13, 1, 2, 3, 26, 5, 2, 0, 4, 7, 1, 6, 1, …
$ antifamily_count      <dbl> 0, 2, 5, 0, 1, 0, 1, 0, 0, 0, 1, 2, 0, 0, 0, 0, …
$ occult_count          <dbl> 0, 3, 2, 0, 0, 0, 3, 1, 0, 0, 0, 1, 1, 1, 0, 0, …
$ language_count        <dbl> 0, 16, 14, 1, 3, 2, 14, 3, 2, 0, 7, 9, 3, 5, 1, …
$ lgbtq_count           <dbl> 1, 5, 1, 0, 2, 0, 10, 1, 0, 0, 3, 3, 0, 0, 0, 0,…
$ violent_count         <dbl> 1, 6, 2, 0, 2, 0, 19, 2, 1, 0, 3, 3, 1, 3, 1, 0,…
$ removed_count         <dbl> 1, 11, 10, 2, 0, 1, 5, 3, 2, 0, 4, 11, 0, 10, 0,…
$ removed_pct           <dbl> 25.000000, 7.432432, 34.482759, 66.666667, 0.000…
$ POPU_LSA              <dbl> 544270, 12054201, 11551941, 1588981, 5633514, 14…
$ POPU_ST               <dbl> 544270, 12284183, 11551941, 2009671, 5618344, 10…
$ CENTLIB               <dbl> 23, 451, 240, 91, 15, 47, 122, 378, 80, 19, 61, …
$ BRANLIB               <dbl> 53, 183, 480, 27, 169, 26, 94, 82, 10, 67, 331, …
$ BKMOB                 <dbl> 2, 29, 58, 1, 18, 2, 8, 6, 12, 4, 17, 6, 1, 4, 7…
$ BKVOL                 <dbl> 2497545, 27790282, 45224425, 4518130, 13954140, …
$ LIBRARIA              <dbl> 189.16, 1445.81, 2725.65, 290.37, 1297.29, 232.1…
$ VISITS                <dbl> 3872783, 47188171, 88255852, 8324986, 33662473, …
$ REGBOR                <dbl> 383126, 5549200, 8767349, 1140637, 3271760, 5747…
bkmob_plot <- bk_df |>
    ggplot(aes(x = BKMOB)) +
    geom_histogram()
bkmob_plot |> class()
[1] "gg"     "ggplot"

Data

tigris_vars <- c("STUSPS10", "NAME10", "ALAND10", "AWATER10", "INTPTLAT10", "INTPTLON10")
main_states <- bk_df |> pull(state)
library(tigris)
state_sf <- states(year = 2010, progress_bar = FALSE) |> 
    select(any_of(tigris_vars)) |>
    filter(STUSPS10 %in% main_states)

state_sf |> st_crs() |> print()
Coordinate Reference System:
  User input: NAD83 
  wkt:
GEOGCRS["NAD83",
    DATUM["North American Datum 1983",
        ELLIPSOID["GRS 1980",6378137,298.257222101,
            LENGTHUNIT["metre",1]]],
    PRIMEM["Greenwich",0,
        ANGLEUNIT["degree",0.0174532925199433]],
    CS[ellipsoidal,2],
        AXIS["latitude",north,
            ORDER[1],
            ANGLEUNIT["degree",0.0174532925199433]],
        AXIS["longitude",east,
            ORDER[2],
            ANGLEUNIT["degree",0.0174532925199433]],
    ID["EPSG",4269]]
bk_sf <- state_sf |>
    left_join(bk_df, join_by(STUSPS10 == state))
bk_sf |> slice_sample(n = 2)
STUSPS10 NAME10 ALAND10 AWATER10 INTPTLAT10 INTPTLON10 division_code division_str political_value_index median_income hs_grad_rate college_grad_rate challenge_count explicit_count antifamily_count occult_count language_count lgbtq_count violent_count removed_count removed_pct POPU_LSA POPU_ST CENTLIB BRANLIB BKMOB BKVOL LIBRARIA VISITS REGBOR geometry
TN Tennessee 106797885992 2355218865 +35.8585639 -086.3493573 6 East South Central -8.7 -2995.5 -3.661958 -4.42373 13 3 0 0 4 2 1 5 38.461539 6270431 6270431 186 103 3 11836785 554.89 21161771 3239146 MULTIPOLYGON (((-88.38051 3…
OR Oregon 248607802255 6191433228 +43.9715225 -120.6226269 9 Pacific 4.0 1274.5 5.538042 1.07627 118 26 1 3 14 10 19 5 4.237288 3657597 3823465 122 94 8 10052777 503.41 25034546 1938999 MULTIPOLYGON (((-124.5523 4…

Aesthetic Mapping

Common aesthetic attributes

  • Position (x, y)
  • Color (or, fill)
  • Other visual features: size and shape (size, shape, linewidth)
bk_sf |>
    ggplot(
        aes(
            x = CENTLIB,
            y = BKMOB,
            color = division_str,
            label = STUSPS10
        )
    ) +
    geom_label(size = 2) +
    labs(
        x = 'Central Library Count',
        y = 'Bookmobile Count'
    )

Tip

labs(): to customize plot labels

Save your plot

plot1 <- bk_sf |>
    ggplot(
        aes(
            x = CENTLIB,
            y = BKMOB,
            color = division_str,
            label = STUSPS10
        )
    ) +
    geom_label(size = 3.5) +
    labs(
        x = 'Central Library Count',
        y = 'Bookmobile Count'
    )
ggsave(
    filename = 'lib-bkmob-scatterplot.png',
    plot = plot1,
    width = 20,
    height = 10,
    dpi = 300
)

geom_sf

bk_sf |> mapview(zcol = 'BKMOB')
bk_sf |>
    ggplot() +
    geom_sf(aes(fill = (BKMOB / POPU_LSA))) +
    geom_sf_label(aes(label = STUSPS10), size = 1.5, fill = "white", alpha = 0.7, color = "black") +
    theme_classic() +
    scale_fill_viridis() +
    labs(
        x = 'Longitude',
        y = 'Latitude',
        fill = "Number of Bookmobiles per Capita"
    )